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Abstract 

Detection  and  modeling  of  self- similarity  and  symmetry  is  important  in 
shape  recognition,  matching,  synthesis,  and  reconstruction.  While  the  detec¬ 
tion  of  rigid  shape  symmetries  is  well  established,  the  study  of  symmetries 
in  non-rigid  shapes  is  a  much  less  researched  problem.  A  particularly  chal¬ 
lenging  setting  is  the  detection  of  symmetries  in  non-rigid  shapes  affected 
by  topological  noise  and  asymmetric  connectivity.  In  this  paper,  we  treat 
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shapes  as  metric  spaces,  with  the  metric  induced  by  diffusion  distances,  and 
define  non-rigid  symmetries  as  self-isometries  with  respect  to  the  diffusion 
metric.  Experimental  results  show  the  advantage  of  the  diffusion  metric  over 
the  previously  proposed  geodesic  metric  for  exploring  intrinsic  symmetries 
of  bendable  shapes  with  possible  topological  irregularities. 


1  Introduction 

Symmetry  or  geometric  self-similarity  plays  a  fundamental  role  in  Nature  and  is 
tightly  tied  to  shape  analysis  and  understanding.  The  exact  definition  of  sym¬ 
metry  heavily  depends  on  what  we  understand  by  the  geometry  of  the  shape. 
Broadly  speaking,  intrinsic  geometry  describes  the  properties  of  the  shape  which 
are  invariant  to  inelastic  deformations,  while  extrinsic  geometry  is  associated 
with  rigid  transformations.  In  the  pattern  recognition  and  computer  vision  lit¬ 
erature,  there  exists  a  significant  number  of  papers  dedicated  to  finding  sym¬ 
metries  in  images  [29],  two-dimensional  [1,  2],  and  three-dimensional  shapes 
[22,  32,  34,  35].  A  wide  spectrum  of  methods  employed  for  this  purpose  includes 
approaches  based  on  dual  spaces  [13],  genetic  algorithms  [17],  moments  calcula¬ 
tion  [10],  pair  matching  [12,  27],  and  local  shape  descriptors  [47]  to  name  just  a 
few.  Traditionally,  symmetries  are  described  as  extrinsic  self- similarity  under  ro¬ 
tations  and  reflections.  While  being  adequate  for  rigid  shapes,  such  a  description 
is  inappropriate  for  non-rigid  ones.  Extrinsic  symmetry  can  be  broken  as  a  re¬ 
sult  of  non-rigid  shape  deformations,  whereas  its  intrinsic  symmetry  is  preserved 
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[37,  38], 

While  finding  shape  representations  that  are  invariant  to  bending  has  been 
intensively  studied  in  recent  years  [5,  14,  15,  18,  19,  28,  30,  31,  40,  42],  very 
little  has  been  done  in  the  research  of  regular  surfaces  and  self-similarity  study 
in  non-rigid  shapes.  Raviv  et  al.  [37]  first  defined  and  proposed  a  computational 
framework  for  exploring  exact  and  approximate  intrinsic  self-similarities  based  on 
invariance  of  geodesic  distances.  They  used  the  numerical  scheme  presented  in 
[5,  6]  for  the  computation  of  minimum-distortion  correspondence  between  metric 
spaces.  In  [34],  Ovsjanikov  et  al.  showed  how  eigenfunctions  of  the  Laplace- 
Beltrami  operator  can  be  used  for  identifying  reflective  intrinsic  symmetries.  They 
transformed  intrinsic  symmetries  to  (approximate)  Euclidean  ones  in  a  feature 
space  created  by  the  Laplace-Beltrami  eigenfunctions.  The  reader  is  also  referred 
to  recent  papers  on  such  spectral  signatures  [3,  24,  39,  40,  42]. 

A  major  drawback  of  the  mentioned  approaches  is  that  topology  and  connec¬ 
tivity  changes  (resulting,  e.g.,  from  noise  and  acquisition  artifacts)  can  have  a 
dramatic  influence  on  the  intrinsic  geometry  defined  by  the  geodesic  distances. 
Since  the  geodesic  is  the  shortest  path,  even  a  small  change  in  the  topology  can 
cause  a  significant  change  in  the  geodesic  metric  by  affecting  the  length  of  many 
geodesics,  thus  the  approach  proposed  in  [37]  is  sensitive  to  topology  changes. 
Ovsjanikov  et  al.  [34]  explain  why  topology  noise  appears  in  high  frequencies, 
hence  its  influence  can  be  reduced  in  their  method.  Yet,  they  can  not  guarantee 
robustness  to  major  topology  changes. 

In  [8,  30]  it  was  shown  that  diffusion  geometry,  arising  from  the  study  of  heat 
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propagation  on  the  surface,  can  gracefully  handle  topological  and  connectivity 
problems .  In  this  paper,  we  introduce  diffusion  symmetries  of  non-rigid  shapes 
which  are  robust  to  topology  changes. 

The  rest  of  this  paper  is  organized  as  follows:  In  Section  2,  we  present  metric 
models  of  shapes.  In  Section  3,  we  define  exact  and  approximate  symmetries,  and 
in  Section  4  elaborate  on  Gromov-Hausdorff  distance.  Section  5  is  dedicated  for 
the  numerical  framework,  followed  by  Section  6  where  we  present  experimental 
results.  Section  7  concludes  the  paper. 

2  Shapes  as  metric  spaces 

We  model  a  shape  as  a  metric  space  {X^dx),  where  X  is  a  two-dimensional  smooth 
compact  connected  and  complete  Riemannian  surface  (possibly  with  boundary) 
embedded  into  M^,  and  d  :X  xX  ^  R+  U  {0}  is  a  metric  measuring  distances  be¬ 
tween  pairs  of  points  on  Z.  In  general,  d  is  independent  of  X  and  can  be  dehned 
quite  arbitrarily,  provided  that  it  satisfies  the  metric  axioms  ,positive  definiteness, 
symmetry,  and  triangle  inequality.  There  are,  however,  a  few  natural  choices  of 
d  describing  different  properties  of  the  shape,  and  giving  rise  to  different  geome¬ 
tries. 
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2.1  Euclidean  geometry 

Considering  X  as  a  subset  of  the  Euclidean  space  E  =  we  can  measure  the 
distances  between  points  x  and  x'  on  X  using  the  restricted  Euclidean  metric, 

d^{x,x)  =  d^\xxx{x,x).  (1) 

We  refer  to  the  properties  described  by  d^  as  the  extrinsic  geometry  of  X.  Such 
properties  are  invariant  under  Euclidean  transformations  ,  translation,  rotation, 
and  reflection  in  E.  However,  nonrigid  deformations  may  change  the  extrinsic 
geometry,  and  thus  the  Euclidean  metric  is  not  suitable  for  modeling  shapes  ex¬ 
hibiting  significant  deformations. 

2.2  Geodesic  geometry 

Geometries  invariant  to  bending,  like  inelastic  deformations  which  do  not  stretch 
or  tear  the  shape,  are  called  intrinsic  and  are  defined  using  properties  of  X  without 
resorting  to  the  ambient  space  E.  The  simplest  example  of  an  intrinsic  metric  is 
the  geodesic  metric,  defined  as  the  length  of  the  shortest  path  on  the  shape, 

dx{x,x)  =  inf  l{y)  (2) 

yer(x,x') 

where  T{x,x')  is  the  set  of  all  admissible  paths  between  the  points  x  and  x'  on  X 
and  £(7)  is  the  length  of  the  path  7  ^  The  geodesic  metric  is  invariant  to  inelastic 

^The  assumption  of  completeness  implies  that  dx  is  realized  by  a  shortest  path  and  thus  replace 
inf  with  min. 


5 


Technion  -  Computer  Science  Department  -  Technical  Report  CIS-2009-20  -  2009 


shape  deformations,  and  as  a  particular  case  is  also  invariant  to  rigid  transfor¬ 
mations.  As  a  consequence,  the  geodesic  metric  is  a  good  choice  for  describing 
non-rigid  shapes  capturing  their  properties  that  do  not  depend  on  deformations. 

2.3  Diffusion  geometry 

Well-known  in  other  fields  for  decades,  diffusion  geometry  has  been  popularized 
in  data  analysis  by  Lafon  and  Coifman  [11,  25].  Informally,  the  diffusion  distances 
have  an  effect  of  averaging  over  all  possible  paths  connecting  two  points,  while 
the  geodesic  distance  is  the  length  of  just  the  minimal  one.  As  a  consequence, 
the  former  is  more  robust  than  the  geodesic  distance  in  cases  where  topological 
changes  are  present  [33,  43]. 

Formally,  the  diffusion  distance,  is  related  to  the  heat  equation, 

M  =  0,  (3) 

governing  the  distribution  of  heat  u  on  the  surface.  Here,  Ax  denotes  the  Laplace- 
Beltrami  operator,  a  generalization  of  the  Laplacian  to  non-Euclidean  domains. 
(Note  that  we  define  Ax  as  a  positive- semidefinite  operator,  hence  its  positive  sign 
in  (3)).  The  fundamental  solution  ht{x,z)  of  the  heat  equation  (3),  also  called  the 
heat  kernel,  is  the  solution  with  a  point  heat  source  at  x  G  A  used  as  the  initial 
condition  at  t  =  0.  The  heat  kernel  ht{x,z)  describes  the  heat  distribution  at  time  t 
at  point  z  e  A. 

For  compact  manifolds,  the  Laplace-Beltrami  operator  has  discrete  eigende- 
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composition  of  the  form 


Ax4>i  =  X4i,  (4) 

where  Ao,Ai,...  are  eigenvalues  and  are  the  corresponding  eigenfune- 

tions.  Using  this  basis,  the  heat  kernel  ean  be  presented  as  [21] 

oo 

ht(x,z)  =  (5) 

i=0 

Note  that  for  a  shape  eonsisteney  of  a  single  eonneeted  eomponent  Aq  =  0,  = 

eonst,  and  A,  >  0. 

The  diffusion  distance  is  the  family  of  metries  of  the  form 

dxMy)  =  (6) 

=  /  \ht{x,z) -ht{y,z)\^dz 

Jx 

oo 

i=i 

parametrized  by  the  seale  t.  Sinee  the  diffusion  distanee  is  derived  from  the 
Laplaee-Beltrami  operator  whieh  is  an  intrinsie  property  of  the  shape,  it  is  an 
intrinsie  metric  and  is  therefore  bending-invariant.  Besides  that,  it  shows  better 
robustness  to  topology  ehanges,  as  we  present  in  the  following. 
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3  Symmetry 


The  metric  model  of  shapes  allows  to  formalize  the  notions  of  invariance,  self¬ 
similarity,  and  symmetry.  From  the  metric  viewpoint,  a  shape  is  invariant  under 
a  transformation  if  the  metric  structure  is  unaffected  by  the  transformation.  More 
formally,  two  shapes  {X,d)  and  (T,  J),  where  d,d  are  the  respective  metrics,  are 
isometric  if  there  exists  a  bijective  map  (p  :  X  Y  (called  isometry)  such  that 
do  [cp  X  cp)  =  d.  In  particular,  the  shape  {X,d)  is  self-similar  or  self-isometric  if 
there  exists  a  permutation  (bijective  map  /  from  X  to  itself)  which  is  an  isometry 
with  respect  to  d,  i.e.,  do  (f  x  f)  =  d. 

It  can  be  easily  shown  that  symmetries  form  a  subgroup  of  the  group  of  permu¬ 
tations  (n(X),o),  where  n(Z)  =  {g  :  and  o  is  the  function  composition 

operator.  We  denote  the  group  of  symmetries  of  X  by 

Sym{X,d)  =  {fe  n(Z)  :do{fxf)=d}.  (7) 

3.1  Approximate  symmetries 

Because  the  true  symmetry  is  an  idealized  rather  than  naturally  occurring  phe¬ 
nomenon,  a  relaxation  of  the  notion  of  symmetry  is  required  [47].  In  the  Euclidean 
case,  approximate  symmetry  can  be  defined  as  a  non-trivial  Euclidean  transforma¬ 
tion  i  satisfying  i{X)  p^X.  In  the  more  general  setting,  the  relaxation  of  the  notion 
of  symmetry  is  possible  by  allowing  d  o  (f  x  f)  ^  d  in  the  definition  of  self¬ 
isometry.  More  formally,  we  say  that  the  shape  {X,d)  is  e-self-isometric  if  there 
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exists  a  permutation  /  G  n(Z)  with  distortion 


dis(/,d)  =  \\do{fxf)-d\\c 


(8) 


=  max  \d{x,x)  —  d{f{x),f{x))\ <  e. 
x,x'ex 


We  denote  the  family  of  all  e-self-isometries  of  {X,d)  by 


Isoe(X,d)  =  {/en(X):dis(g,/)<e}. 


(9) 


In  particular,  for  e  =  0  we  get  the  symmetry  group  Isoo(X,<i)  =  Iso(Z,  J).  Note 
that  unlike  Iso(X,<i),  IsOe(X,<i)  is  not  a  group  for  e  >  0,  since  it  does  not  satisfy 
the  closure  property:  a  composition  of  two  e-self-isometries  is  in  general  a  2e- 
self-isometry. 


3.2  Local  asymmetry 

In  many  cases,  the  shape  asymmetry  is  as  important  as  symmetry,  as  it  indicates 
some  abnormality  of  a  shape.  Given  an  e-symmetric  shape  and  one  of  its  approx¬ 
imate  symmetries  /  E  IsOe(Z,  J),  we  can  compute  the  contribution  of  each  point 
to  the  shape  asymmetry  as  the  local  shape  asymmetry, 

=  max\dx{x,x’)-dx{f{x)J{x’))\,  (10) 

X 

quantifying  the  distortion  of  /  at  a  point  x.  Points  with  large  local  asymmetry  are 
responsible  for  symmetry  breaking.  Global  shape  asymmetry  can  be  calculated  as 
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the  supremum  of  all  local  ones,  which  is  the  same  as  the  distortion  dis(/,<i). 

3.3  Symmetry  space 

Though  we  cannot  use  a  group  structure  to  represent  approximate  symmetries, 
we  can  still  explore  the  space  of  permutation  n(X),  where  each  function  g  has 
its  distortion  dis(g,<i).  If  the  distortion  is  a  smooth  function  then  choosing  local 
minima  as  symmetry  candidates  is  the  logical  thing  to  do.  For  smooth  mappings 
this  is  the  case,  which  means  that  local  minima  of  distortions  are  good  candidates 
for  signihcant  approximate  symmetry.  See  supporting  proof  in  Appendix  A  for 
the  continuity  of  the  distortion  near  true  symmetries. 

The  space  of  functions  n(X)  can  also  he  endowed  with  some  standard  metric 
that  measures  the  distance  between  two  permutations  /  and  g  on  A,  defined  as 

^n(x)(C,?)  =maxcf(/(x),g(x))  =  d{f{X),g{X)),  (11) 

which,  in  turn,  depends  on  the  choice  of  the  metric  d.  We  refer  to  the  set 

^n(x) =  {/  ^  n(A)  :  dii(^x)i8if)  <  (12) 

as  the  closed  metric  ball  (intrinsic  or  extrinsic,  according  to  the  choice  of  the 
metric  in  the  dehnition  of  ^(n(x))  of  radius  r  centered  at  g.  (We  will  omit  r  referring 
to  a  ball  of  some  unspecified  radius). 

Since  true  symmetries  have  zero  distortion  (i.e.  dis(g,  J)  =  0),  they  are  the 
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global  minimizers  of  the  distortion  function  dis(g,(i).  Moreover,  they  are  also 
local  minimizers  of  the  distortion,  in  the  sense  that  for  every  approximate  sym¬ 
metry  g,  there  exists  a  sufficiently  small  neighborhood  5n(x)(g)>  such  that  any 
/  ^  ^n(x)(.?)’  dis(/)  >  dis(g).  While  trivial  for  true  symmetries,  this  observa¬ 
tion  motivates  our  definition  of  approximate  symmetries  as  local  minima  of  the 
distortion. 

We  therefore  define  approximate  symmetries  as 

Syme(X,cf)  =  {g  e  Isoe(Z,cf)  :  dis(g,cf)  <  dis(/,d)V/  G  5n(jf)(g)}.  (13) 

Numerically,  we  observed  that  for  a  small  e,  the  number  of  local  minima  is  similar 
to  the  number  of  true  symmetries  of  a  symmetric  shape. 

4  Gromov-Hausdorff  distance 

An  early  attempt  to  compare  shapes  as  metric  spaces  was  done  by  Elad  and  Kim- 
mel  [14]  who  proposed  to  embed  pairwise  geodesic  distances  into  a  third  fixed 
metric  space  {Z,dz)-  Usually,  a  low  dimensional  Euclidean  space  is  a  good  candi¬ 
date  from  the  computational  point  of  view.  Such  a  minimum  distortion  embedding 
was  referred  to  as  canonical  form.  It  can  be  explicitly  calculated  by  minimizing 

min  mdix\dx{x,x)  —  dz{^{x),^{x)\.  (14) 

(p:X-^Zx,x'eX 
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Assuming  the  embedding  distortion  is  low  enough,  the  comparison  of  {X,dx) 
and  {Y^dy)  can  be  reduced  to  rigid  matching  of  canonical  forms,  e.g.  using  the 
Hausdorff  distance. 


i//(y))  =  max{  max  min  dz{x,y),  max  min  dz{x,y)}.  (15) 

yeY{Y)xe(p{x)  xe(p{x)ye'Hf{Y) 


A  general  theoretical  framework  for  metric  spaces  comparison  was  proposed  by 
Gromov  [9,  18],  and  introduced  to  shape  analysis  by  Memoli  ans  Sapiro  [31]. 
Since  it  is  generally  impossible  to  select  a  common  metric  space  {Z,dz)  accom¬ 
modating  all  shapes,  Gromov  suggested  to  optimize  also  for  Z,  resulting  in  the 
following  distance,  referred  to  as  the  Gromov-Hausdorff  distance. 


dGH{X,Y)=  inf  4((p(A),v/(y)),  (16) 

(p:X~*Z 

Z 

where  y/  and  <p  are  isometric  mappings.  The  Gromov-Hausdorff  distance  can 
be  alternatively  expressed  in  terms  of  correspondences  between  the  two  metric 
spaces  [9]  as 


dGHiX,Y)  =  ^mfdisiC), 


(17) 
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where  C  C  Z  x  F  is  a  correspondence  between  the  spaces,  and  dis(C)  is  its  distor¬ 
tion 


dis(C)  =  sup  \dx{x,x)-dY{y,y)\.  (18) 

{x,y),{x’,y’)eC 

Note  that  for  C  to  be  a  valid  correspondence,  for  each  x  G  X  there  must  exist  at 
least  one  y  G  F  such  that  (x,  y)  G  C,  and  vice  versa,  for  each  y  G  F  there  must  exist 
at  least  one  x  G  X  such  that  (x,y)  G  C. 

In  [5],  it  was  shown  how  (17)  can  be  efficiently  approximated  using  a  convex 
optimization  framework.  In  this  paper,  we  adopt  this  scheme  for  the  computation 
of  self-isometries. 


5  Numerical  framework 

For  practical  computation  of  symmetries,  the  surface  X  has  to  be  discretized  and 
sampled  at  N  points,  constituting  an  r-covering  (i.e.,  X  =  [j^=iBx{x„,r),  where 
Bx  denotes  a  closed  ball  on  X  induced  by  the  metric  dx).  We  denote  this  sampling 
by  Xr  =  {xi,  ...,XAf}  C  X.  A  good  sampling  strategy  can  be  achieved  using  the  (2 
optimal)  farthest  point  sampling  algorithm  [20],  which  guarantees  that  X^  is  also 
r-separated,  i.e.  dx,t{xi,Xj)  >  r  for  any  i  f  j.  In  the  following,  we  assume  thatX 
is  approximated  by  a  triangular  mesh  X  built  upon  the  vertices  X^. 
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5.1  Approximation  of  distances 

In  order  to  approximate  the  geodesic  distances,  we  use  the  fast  marching  method 
[23],  which  produces  a  first-order  approximation  for  the  geodesic  distances  be¬ 
tween  Xr  on  Z.  The  main  idea  of  fast  marching  is  to  simulate  wavefront  propa¬ 
gation  on  a  triangular  mesh,  associating  the  time  of  arrival  of  the  front  with  the 
distance  traveled,  assuming  constant  propagation  speed.  A  linear  time  is  required 
to  obtain  the  geodesic  distance  between  two  points  [44,  45,  46]. 

The  diffusion  metric  is  computed  according  to  (6).  For  this  purpose,  we  first 
compute  the  discrete  approximation  of  the  Laplace-Beltrami  operator  on  the  mesh. 
A  Finite  Elements  Method  (FEM)  approach  is  possible  [39],  but  in  this  work  we 
used  the  following  generic  form 


(%/).'  = (19) 

at  j 

where  /  :  A  — )>  R  is  a  scalar  function  defined  on  the  mesh  X  and  represented  as  a 
vector  of  function  values  at  the  vertices  of  the  mesh,  Wij  are  weights,  and  ai  are 
normalization  coefficients.  In  matrix  notation.  Equation  (19)  can  be  written  as 

A^/  =  A-T/,  (20) 


where  A  =  diag(a/)  and  L  =  diag  iXi^iWu)  -  {wij). 

There  exists  several  discretizations  of  the  Laplace-Beltrami  operator  [4,  48], 
for  which  A  and  W  are  defined.  Here  we  adopt  the  cotangent  weight  scheme 
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[36],  in  which  at  is  set  to  be  half  of  the  area  of  the  eircumeentrie  dual  of  vertex  i, 
and  Wij  =  eota,j  -l-eotj8/;,  where  aij  and  are  the  angles  opposite  to  the  edge 
between  the  i’th  and  /th  vertiees,  if  they  are  eonneeted,  and  zero  otherwise. 

By  solving  the  generalized  eigendecomposition  problem  [26] 

A(t>  =  XL(t),  (21) 

the  k  smallest  eigenvalues  Ao, and  eorresponding  eigenfunetions  •••5  : 

X  — ^  M  of  the  diseretized  Laplaee-Beltrami  operator  are  computed. 

5.2  Coarse  matching 

Diserete  permutations  on  can  be  represented  as  X-tuples  of  the  form  n  = 
e  {1,...,X}^.  Finding  all  permutations  with  a  distortion  lower  than 
e  requires  eomputation  of  the  distortion  of  id'(N\)  possible  mappings.  However, 
the  seareh  spaee  can  be  greatly  reduced  by  ruling  out  mappings  that  cannot  have 
low  distortion.  We  observe  that  for  a  good  candidate  of  an  approximate  symmetry, 
the  intrinsie  properties  of  the  surfaee,  such  as  local  intrinsic  geometry  around  xi 
should  be  similar  to  that  around  In  order  to  quantify  this  behavior,  for  eaeh 
Xi  G  Xr  we  eompute  the  histogram  hi  =  hist{{dij  :  dij  <  p})  of  the  approximate 
distanees  (dij  =  dx{xi,Xj))  in  a  p-ball  centered  at  x,  [16,  28,  37]. 

In  order  to  reduee  the  number  of  possible  permutations,  we  ehose  a  small 
number  of  points  using  farthest  point  sampling  strategy  on  diffusion  distances  to 
represent  the  symmetry.  Experiments  showed  that  resampling  the  surfaee  starting 
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from  the  second  sampled  point  produces  a  better  sampling  because  the  location  of 
the  first  point  is  arbitrary  and  its  support  may  occlude  a  better  one.  Each  point  has 
many  possible  matches  and  clustering  of  the  candidates  is  performed.  Two  points 
are  in  the  same  cluster  if  and  only  if  there  exists  a  chain  of  points  on  the  surface 
with  similar  distance  histograms.  Earth  mover’s  distance  (EMD)  [41]  can  be  used 
as  the  distance  between  histograms. 

Eor  a  given  point  x,  we  compute  its  histogram  hx  and  find  all  points  y  with  a 
similar  histogram,  denoted  by 


H5W  =  {y  ■  d{hx,hy)  <  5},  (22) 

where  d  is  the  distance  between  the  histograms,  and  5  is  calculated  a  priori  ac¬ 
cording  to  the  histogram  distance  between  neighboring  points.  We  perform  e 
separation  of  the  set  (x)  to  subsets  (x)  according  to 

y,  z  e  Hg  (x)  iff  y,  z  e  H5  (x)  and  3  (y^jf  G  Us  (x) 

s.t.  yi  =  y,  y™  =  z,  VI  <  <  m  -  1  (23) 

dx{yk,yk+i)  <  £■ 

Once  possible  sets  are  constructed,  we  search  for  a  global  correspondence 
between  candidates  of  each  set.  Since  the  number  of  sets  is  small,  the  search 
becomes  feasible.  Solving  (23)  is  done  in  similar  spirit  to  union-set  data  structure 
methods,  as  follows: 
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1.  Choose  one  point  y  at  random. 

2.  Find  the  farthest  point  q  from  it. 

3.  Resample  the  mesh  \x%mg  farthest  point  sampling  starting  with  q. 

4.  Repeat  the  following  stages  for  eaeh  point  sample  x. 

(a)  Calculate  Hist^  (jc) . 

(b)  i=l.  V/  Si  =  {} 

(c)  Choose  a  point  z  G  Hist^  (x)  \  Ufe=  i  'S'fe  • 

(d)  Set  Si  =  {z}. 

(e)  For  each  point  q  E  Sf. 

Si  =  e  H5(x)  s.t.  d{q,g)  <  e} 

(f)  repeat  step  (e)  until  convergence. 

(g)  if  ULi  Sk  7^  H5(x):  i=i+l,  goto  step  (c) 

(h)  =  argmm^{d{hx,hq)  s.t.  <?  G  Sj. 

Each  point  x  has  {mf }  set  as  possible  matches. 

In  many  practical  problems,  exhaustive  search  is  possible.  In  cases  where 
search  complexity  is  still  prohibitive,  we  use  a  branch  and  bound  algorithm  similar 
to  [16].  While  complexity  remains  the  same,  the  search  space  becomes  extremely 
small  in  practice.  As  presented  in  [37],  we  apply  a  branch  and  bound  procedure, 
for  M  points,  as  follows. 
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1 .  Given  a  correspondence  ofk—l  feature  points  l)i— Kk- 1 ) , 
we  would  like  to  establish  k  ^  TZk- 

2.  Prune:  For  each  potential  correspondence  Ttk  G  Q,  evaluate  max/=i^ ^ \dik  — 
dTti,Ttk\-  If  obtained  distortion  is  larger  than  e,  discard  the  potential  cor¬ 
respondence. 

3.  Branch:  For  each  remaining  n^,  recursively  invoke  Step  1  with  (1,  ...,k)  i— > 

4.  Bound:  If  k  =  M,  compute  the  distortion  dis(;r,<i).  If  dis(;r,<i)  <  e,  mark 
n  as  an  e- symmetry. 

5.3  Generalized  MDS 

Once  a  coarse  match  is  found,  it  is  used  as  an  initialization  of  the  second  rehne- 
ment  stage.  Solving  (17)  requires  to  hnd  a  correspondence  which  distorts  the 
metric  the  least.  Following  the  discretization  in  [5],  for  N  mesh  points  Xi  on  X 
we  search  for  x'  locations,  (not  necessarily  coincide  with  the  mesh  vertices),  with 
similar  diffusion  distances  between  them,  minimizing 


I  /s  "^  /  /  /  \  1 2 

min  max  \di;  —  dxiXi.xAi  .  (24) 

We  use  the  convex  optimization  scheme  presented  in  [5],  referred  to  as  generalized 
multidimensional  scaling  (GMDS),  to  solve  (24). 
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6  Experimental  results 

In  order  to  show  the  advantages  using  diffusion  distances,  we  performed  sev¬ 
eral  experiments  on  meshes  taken  from  the  TOSCA  dataset  [7],  The  first  coarse 
matching  exploration  stage  took  a  couple  of  seconds  per  surface,  and  the  second 
stage  using  GMDS  rehnement  procedure  took  less  than  a  minute  for  approxi¬ 
mately  2000  points  on  a  Pentium  Core  2  Duo  3.0  GHz  CPU. 

Figure  1  depicts  the  histograms  of  geodesic  and  diffusion  distances  of  two 
matching  points.  Without  a  topology  change,  all  histograms  are  almost  identical, 
while  with  a  topology  change  only  the  diffusion  histograms  preserve  their  simi¬ 
larity. 


Diffusion  Histograms  Geodesic  Histograms 

With  Topology  changes  With  Topology  changes 


Without  Topology  changes  Without  Topology  changes 


Figure  1 :  Geodesic  and  diffusion  histograms  of  two  matching  sample  points  with 
and  without  topology  changes.  Black  circles  mark  points  where  topology  changes 
were  introduced,  and  over  the  major  non-matching  parts  of  the  histograms. 


In  Figure  2  we  present  geodesic  and  diffusion  maps  from  a  point  source.  In 
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contrast  to  geodesic  distances,  diffusion  distanees  are  largely  unaffeeted  by  topol¬ 
ogy  ehanges. 


Without  Topology  change  With  Topology  change 


Geodesics 


Diffusion 


Figure  2:  Color  eoded  distanees  from  the  foot  (marked  with  an  arrow)  using  two 
different  metrics  under  topology  ehanges.  In  contrast  to  geodesic  distances,  diffu¬ 
sion  distances  are  hardly  affeeted  by  topology  ehanges.  Blaek  cireles  mark  points 
where  topology  changes  were  introduced. 


Figure  3  shows  the  best  five  eoarse  symmetries  of  a  human  body  undergoing 
several  topology  ehanges.  The  reflective  symmetry  (E)  is  found  as  well  as  ad¬ 
ditional  twists  (B,C,D).  In  a  small  number  of  points  false  negatives  appear,  and 
these  ean  be  rejeeted  using  either  the  GMDS  proeedure,  adding  descriptors,  or 
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increasing  the  sample  size. 


0 


0.015351 


0.010598 


0.015654 


0.016473 


Figure  3:  Coarse  symmetries  of  a  human  body.  Similar  parts  appear  in  similar 
eolors,  and  the  total  distortion  is  indieated  as  a  numerieal  value.  The  distortion 
less  identity  symmetry  is  (A).  The  full  refleetive  symmetry  is  (E)  while  (B),  (C), 
and  (D)  represent  possible  symmetries  with  partial  reflection.  Black  circles  mark 
points  where  topology  changes  were  introduced. 


The  priee  paid  to  eope  with  topology  ehanges  is  diseriminatively.  One  can 
think  of  diffusion  distances  as  an  averaging  process.  As  sueh  it  diseriminates  less 
the  geometrie  features.  From  Figure  4  it  is  evident  that  the  diffusion  histograms 
of  the  palms,  the  legs  and  the  faee,  are  similar.  Measuring  the  distance  between 
the  diffusion  histogram  descriptors  alone  ean  produce  false  positive  matehings. 
Indeed,  the  histograms  are  not  identieal,  but  to  cope  with  reasonable  deformations 
we  need  to  increase  the  thresholds,  which  can  produce  false  positives.  For  exam- 
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pie,  if  point  4  and  5  should  be  considered  a  match  then  point  3  could  be  a  candidate 
as  well.  Such  false-positives  can  be  removed  using  the  branch-and-bound  proce¬ 
dure,  as  described  earlier,  or  even  by  exhaustive  search  on  all  5 !  permutations,  by 
measuring  the  diffusion  distances  between  those  five  points. 


Histogram  similarity 


Figure  4:  Diffusion  histograms  of  different  points  on  the  body  are  not  identical 
but  similar  due  to  the  diffusion  averaging  process.  A  match  between  points  4  and 
5  can  inflict  a  false-positive  match  with  point  3  as  well.  Black  circles  mark  points 
where  topology  changes  were  introduced. 


Figure  5  depicts  symmetries  constructed  by  using  GMDS.  We  used  L2  instead 
of  Loo  in  the  process,  which  provided  better  numerical  results.  Each  blue  line  rep¬ 
resents  a  possible  match  from  a  set  of  128  matches.  We  show  several  matches 
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on  semi-transparent  human  bodies.  We  calculated  self-matching  on  three  differ¬ 
ent  meshes  with  and  without  topology  changes,  of  different  strengths  (A-C),  us¬ 
ing  geodesic  and  diffusion  distances.  Four  experiments  were  performed  for  each 
mesh,  and  the  quality  of  the  correspondence  with  and  without  a  topology  change, 
for  geodesic  and  diffusion  distances,  was  measured  in  terms  of  the  geodesic  dis¬ 
tances  from  the  ground  truth.  For  ground  truth  x,-  points  and  Xi  calculated  positions 
we  quantified  the  correspondence  quality  by 


d<^{x,x) 


Hjdxjxi^Xi) 
n  ■  Diam(X)  ’ 


(25) 


where  Diam(Z)  is  the  diameter  of  X  and  n  is  the  number  of  points  used  in  the 
optimization  stage,  which  produces  a  scaling  and  sampling  invariant  measure. 

Correspondence  quality  is  summarized  in  Table  1 .  Without  topology  changes 
both  diffusion  distances  and  geodesic  distances  produce  good  results,  but  with 
the  presence  of  topology  changes  the  performance  of  geodesic  distances  is  highly 
degrades  while  diffusion  distances  keep  their  accuracy.  For  example,  in  mesh  (C) 
the  topology  change  caused  a  decrease  in  correspondence  accuracy  by  33%  using 
diffusion  distances,  and  by  more  than  4  times  when  using  geodesic  distances. 

Figure  6  depicts  the  local  shape  asymmetry  of  a  human  body  with  local  asym¬ 
metry  introduced  by  deforming  one  of  the  palms  followed  by  adding  a  topology 
change,  marked  with  a  black  circle.  Even  though  the  body  does  not  exhibit  non¬ 
trivial  extrinsic  symmetries,  we  were  still  able  to  find  its  local  intrinsic  asymmetry, 
marked  with  an  arrow.  High  asymmetry  values  are  shaded  in  blue.  We  extracted 
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same  topology 

different  topology 

Geodesic 

Diffusion 

Geodesic 

Diffusion 

A 

0.0091 

0.0094 

0.0156  (71%) 

0.0126  (34%) 

B 

0.0187 

0.0125 

0.0243  (30%) 

0.0125  (11%) 

C 

0.0160 

0.0222 

0.0901  (463%) 

0.0222  (33%) 

Table  1:  Correpondense  quality  measured  by  d^{x,x)  (25)  on  three  different 
meshes  with  and  without  a  topology  change,  using  geodesic  and  diffusion  dis¬ 
tances.  We  added  the  percentage  decrease  in  accuracy,  due  to  topology  changes, 
in  parentheses.  The  meshes  can  be  seen  in  Figure  5. 


A 


B  C 


Figure  5:  GMDS  self-embeddings.  Matching  points  are  connected  with  a  blue 
line  on  partial  transparent  human  bodies.  (A)  poses  small  topology  changes,  (B) 
mediocre  ,  and  (C)  extreme  ones.  Black  circles  are  drawn  over  connectors  which 
caused  topology  changes. 
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a  non-trivial  intrinsic  symmetry  that  appeared  as  an  intrinsic  reflection,  and  then 
used  that  mapping  to  detect  the  abnormality  at  the  palm.  The  method  had  diffi¬ 
culties  locating  asymmetric  parts  using  geodesic  distances,  yet  succeeded  using 
diffusion  distances. 


Figure  6:  Calculating  local  asymmetry  values  for  a  possible  reflective  symmetry, 
using  geodesic  distances  (right)  and  diffusion  distances  (left).  The  right  palm  was 
stretched  and  a  topology  change  was  inflicted  on  the  legs.  A  black  circle  repre¬ 
sents  the  inflicted  topology  change,  and  an  arrow  represents  the  high  asymmetry 
area  of  the  object.  Blue  color  represents  local  asymmetry.  The  asymmetry  in  the 
palm  was  detected  using  diffusion  distances  (left)  but  could  not  be  isolated  using 
geodesic  distances  (right). 


Figure  7  presents  an  extrinsic  reflective  symmetry  calculated  from  pairs  of 
matching  points.  We  assumed  the  shape  was  extrinsically  symmetric,  and  calcu¬ 
lated  the  reflective  plane  using  Principle  Component  Analysis  (PCA)  on  the  cen¬ 
ters  of  the  connecting  lines  between  matched  points.  Adding  a  topology  change 
rotated  the  reflective  plane  by  almost  10  degrees  when  geodesic  distances  were  in 
use,  but  only  by  one  degree  while  using  diffusion  distances. 
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Figure  7 :  Extrinsic  reflective  symmetry  calculated  from  matching  pairs  of  points 
using  geodesic  (left)  and  diffusion  (right)  distances.  A  topology  noise  is  marked 
in  a  black  circle,  and  matching  points  by  a  red  line.  Compared  to  an  accurate 
reflective  plane,  the  topological  noise  rotated  the  plane  by  9.8  degrees  while  using 
geodesic  distances,  and  by  1.2  degrees  for  diffusion  distances. 
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For  a  given  metric  d,  which  is  t  time  depended,  and  a  given  e,  we  can  sample 
the  space  of  all  possible  symmetries  Sym£(Z,  J).  We  measure  the  distance  be¬ 
tween  mappings  according  to  (11),  and  embed  this  abstract  space  into  using 
multidimensional  scaling  (MDS)  method. 

Since  the  space  is  obviously  non-Euclidean,  there  are  distortions  in  the  pro¬ 
cess,  yet,  the  visual  results  are  informative.  In  Figure  8  we  continued  sampling  the 
symmetries  of  a  human  body  from  Figure  3.  We  computed  155  potential  candi¬ 
dates,  calculated  their  distortion,  embedded  them  in  using  classical  MDS,  and 
interpolated  their  distortion  on  a  plane.  We  see  that  the  symmetry  space  itself  has 
one  approximate  reflective  symmetry,  indicating  that  a  reflective  symmetry  is  the 
only  true  symmetry  of  the  shape.  In  addition,  in  Figure  9  we  present  the  influence 
of  £  on  the  symmetry  space. 


27 


Technion  -  Computer  Science  Department  -  Technical  Report  CIS-2009-20  -  2009 


« 


Figure  8:  The  symmetry  spaee  of  a  the  human  body  from  Figure  3  embedded  in 
Red  represent  low  distortion  embeddings  (symmetries).  The  red  line  shows 
the  refleetive  symmetry  of  the  space  itself.  We  measured  the  distance  between 
functions  according  to  (11). 


£  =  0.1  £  =  0.05 


£  =  0.0125 


I' 


_ fi. 


Figure  9:  £  influence  on  the  symmetry  space  of  the  human  body  from  Figure  5 
embedded  in  R^.  As  e  grows,  more  symmetries  appear,  but  the  symmetry  space 
remains  symmetric. 


28 


Technion  -  Computer  Science  Department  -  Technical  Report  CIS-2009-20  -  2009 


7  Conclusions 

We  presented  a  method  to  extract  approximate  intrinsic  symmetries  of  bendable 
surfaces  which  are  topology  aware,  using  diffusion  distances.  We  showed  how 
local  asymmetry  values  can  be  used  for  detecting  shape  abnormalities,  and  pre¬ 
sented  a  visualization  of  the  entire  symmetry  space. 
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Appendix  A 

Proposition  1 

Let  g  e  n(A)  be  a  continuous  function  of  the  metric  space  {X,dx),  then  dis(g,  J) 
is  a  continuous  function  as  well. 

Proof 

Since  g  is  continuous,  Ve  35  x,x  e  X  such  that  dx{x,x)  <  5 
dx{g{x):g{x)  <  e. 

In  order  to  prove  that  dis  is  continuous  it  is  sufficient  to  show  that 
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Ve  35  x,x^y,y  G  X 


dx  (x,  i)  +  dx  (y,  y)  <  5  (26) 

I  max|Jx(x,y)-rfjf(g(x),^(y))|-max|6?jf(x,y)-cf(g(x),g(y))|  |  <  e. 
x,y  x,y 


But  if  g  is  e-symmetry, 


max  \dx (x, y)  -  dx (g(x) ,g(y))j-  max  jdx (x, y)  -  dx {g{x) ,  g(y) )  | 

x,y  x,y 

<  maxj  \dxix,y)-dxigix),g{y)\-\dx{x,y)-dxigix),g{y))\  \ 
x,y,x,y 

<  maxjdx(x,y)  -£((x,y)|  +  |£(x(g(x),g(y))  -  Jx(g(x),g(y))| 
x,y,x,y 

<25  +  e  +  2e.  (27) 


Then,  for  a  small  e,  and  a  locally  small  e 


5 


e  -  e  -  2e 
2 


(28) 


is  required  for  (26)  to  hold.  Note,  that  e  can  he  infinitely  small  only  near  true 
symmetries.  Obviously,  it  does  not  mean  that  near  e-symmetry  there  exists  an 
^e-symmetry. 
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